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The  photoelectron  spectrum  in  the  ultra-relativistic  limit  of  tunneling  ionization  is  strongly  af¬ 
fected  by  wave-particle  resonance  and  finite  spot-size  effects,  in  contradistinction  with  the  usual 
assumptions  of  strong  field  physics.  Near  term  laser  facilities  will  access  a  regime  where  ionized 
electrons  are  abruptly  accelerated  in  the  laser  propagation  direction,  such  that  they  stay  in  phase 
with  the  laser  fields  through  a  substantial  portion  of  the  confocal  region.  The  final  momentum  of 
the  electron  depends  significantly  on  where  in  the  confocal  region  it  originated.  By  manipulating 
the  target  and  collection  geometry,  it  is  possible  to  obtain  low  emittance,  low  energy  spread,  gi¬ 
gavolt  photoelectrons.  Radiation  reaction  effects  play  a  negligible  role  in  near  term  scenarios,  but 
become  interesting  in  the  multi-exawatt  regime.  A  significant  advance  in  numerical  particle  tracking 
is  introduced. 


I.  INTRODUCTION 

As  laser  technology  continues  to  advance,  new  regimes 
of  laser- matter  interaction  appear.  High-intensity  lasers 
are  associated  with  the  threshold  for  tunneling  ionization 
of  outer  shell  electrons,  about  1014  W/cm2.  Ultra-intense 
lasers  are  associated  with  the  threshold  for  relativistic 
electron  dynamics,  about  1018  W/cm2.  Extreme  light 
facilities  will  produce  intensities  of  1022  W/cm2.  The 
physics  defining  this  regime  is  yet  to  be  determined.  This 
paper  points  out  a  significant  effect  that  is  likely  to  be 
observed:  the  transition  of  photoelectron  spectra  from 
a  ponderomotive  to  a  wave-particle  resonance  regime. 
Fig-  m  illustrates  the  scenario  under  consideration. 

A  fundamental  high  field  process  is  multi-photon  ion¬ 
ization.  The  photoelectron  spectrum  typically  ranges 
from  the  few  eV  range,  as  in  above  threshold  ionization 
(ati)  pa  ei,  to  the  few  MeV  range,  as  in  laser  ioniza¬ 
tion  and  ponderomotive  acceleration  (LIPA)  [3j  4|.  The 
high  order  limit  of  multi-photon  ionization  is  tunneling, 
in  which  the  probability  current  of  the  outgoing  electron 
depends  only  on  the  instantaneous  field  EHE].  Although 
this  leads  to  an  electron  distribution  that  is  initially  con¬ 
fined  to  a  small  region  of  phase,  during  the  subsequent 
motion  electrons  are  exposed  to  a  large  number  of  opti- 
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FIG.  1:  Configuration  of  photoelectron  generation.  The  red 
area  represents  the  confocal  region  of  an  extreme  light  laser 
pulse  propagating  from  left  to  right.  Electrons  ionized  from 
low  charge  states  (black  orbit)  are  spawned  outside  the  ex¬ 
treme  held  region  and  are  ejected  ponderomotively.  Electrons 
ionized  from  high  charge  states  (blue  orbit)  are  spawned  in  an 
extreme  held  and  can  become  resonant  with  the  optical  wave. 
Resonant  electrons  are  forward  directed  and  may  obtain  GeV 
energies. 


cal  cycles  before  leaving  the  confocal  region.  A  new  type 
of  photoelectron  spectrum  is  obtained  when  the  free  elec¬ 
tron  stays  confined  to  a  small  region  of  phase  at  all  points 
in  the  confocal  region,  i.e.,  when  it  is  in  phase  resonance 
with  the  optical  wave.  In  this  paper  the  laser  pulse  pa¬ 
rameters  needed  to  access  this  new  regime  are  clearly 
identified,  and  it  is  shown  that  near-term  multi-petawatt 
laser  facilities  [9]  are  suitable  for  an  experimental  demon¬ 
stration. 

This  work  makes  contact  with  both  the  strong  field 
physics  of  atoms  mg,  and  the  acceleration  of  electrons 
by  electromagnetic  fields  in  free  space  fT0MT5].  Of  the 
many  articles  pertaining  to  strong  field  physics,  those 
concerning  the  photoelectron  spectra  produced  in  the 
LIPA  scheme  13  S]  are  most  closely  related  to  this 
work.  The  process  described  herein  can  be  thought  of 
as  the  extreme-light  limit  of  LIPA,  hereafter  abbreviated 
xLIPA.  In  LIPA,  electrons  originate  from  a  tenuous  gas 
of  moderately  heavy  elements.  Upon  focusing  an  ultra¬ 
intense  laser  pulse  into  the  gas,  electrons  are  tunnel  ion¬ 
ized  and  accelerated.  When  the  ionization  potential  is 
matched  to  the  laser  power  and  focusing,  such  that  the 
electrons  are  ionized  only  when  they  are  exposed  to  the 
peak  laser  intensity,  free  electrons  are  spawned  in  a  large 
ponderomotive  potential,  and  gain  high  energy.  In  the 
case  of  xLIPA,  the  ionization  potential  has  to  be  matched 
in  a  similar  way,  but  the  acceleration  is  no  longer  pon¬ 
deromotive.  Instead,  electrons  stripped  from  atoms  with 
certain  favorable  initial  positions  are  brought  abruptly  to 
the  speed  of  light  and  gain  further  energy  as  they  stay 
in  nearly  the  same  phase  of  the  laser  field  for  extended 
periods  of  time. 

Of  the  many  articles  pertaining  to  free  space  accel¬ 
eration,  those  concerning  the  “Capture  and  Accelerate 
Scenario”  (CAS)  [12J  [14]  are  most  closely  related  to  this 
work.  In  this  scheme,  electrons  are  externally  injected 
into  a  laser  focus,  and  are  abruptly  accelerated  to  the 
speed  of  light.  CAS  differs  from  xLIPA  primarily  in  the 
nature  of  the  initial  conditions.  Namely,  in  xLIPA  the 
intial  conditions  are  tied  to  the  strong  field  physics  of 
atoms,  whereas  in  CAS  they  are  engineered  using  a  con¬ 
ventional  accelerator.  In  practical  terms,  xLIPA  accesses 
the  CAS  regime  without  the  need  for  a  conventional  ac- 
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FIG.  2:  Maximum  energy  in  the  photoelectron  spectrum  for  a 
0.8  fim  wavelength,  20  fs  pulse,  focused  to  a  5  fim.  spot  (pulse 
metrics  are  referenced  to  1  /e  of  the  field) . 

celerator. 


II.  OVERVIEW  OF  RESULTS 

The  primary  result  of  this  paper  in  shown  in  Fig.  [2j 
which  plots  the  maximum  photoelectron  energy  as  a  func¬ 
tion  of  the  normalized  vector  potential,  a  =  eA/mc2. 
The  parameter  a  characterizes  the  importance  of  rela¬ 
tivistic  effects,  i.e.,  when  a  <  1  the  free  electron  dy¬ 
namics  are  weakly-relativistic,  while  when  a  1  they 
are  ultra-relativistic.  Each  value  of  a  is  paired  with  a 
particular  ionization  potential,  [fion,  shown  on  the  upper 
horizontal  axis.  The  U{on  values  are  chosen  so  that  the 
corresponding  a  is  twice  the  threshold  value  for  tunneling 
ionization.  The  plot  shows  simulated  data  with  and  with¬ 
out  radiation  reaction  (RR),  along  with  three  theoretical 
curves  for  comparison.  The  numerical  methods  used  for 
the  simulation  are  in  the  appendix,  and  the  derivation  of 
the  theoretical  curves  is  discussed  below. 

Fig.  §  illustrates  several  important  points.  First, 
two  regimes  of  photoelectron  generation  are  clearly  sug¬ 
gested,  one  for  a  <  100  and  one  for  a  >  100.  The  sim¬ 
ulation  results  for  a  <  10  are  well  matched  by  the  theo¬ 
retical  curve  marked  “ponderomotive,”  which  is  just  the 
prediction  of  LIPA.  The  phase  resonance  model,  which 
is  derived  below,  gives  two  different  curves  depending  on 
whether  axial  field  components  are  kept  or  dropped  in 
the  analysis.  While  both  of  these  curves  fit  the  data  for 
a  >  100  better  than  the  LIPA  model,  the  one  that  ac¬ 
counts  for  the  axial  field  components  is  more  accurate.  In 
simulations  where  the  axial  field  is  artificially  suppressed, 
the  maximum  energy  gain  is  indeed  reduced  as  the  theo¬ 
retical  curves  indicate. 

Perhaps  the  most  interesting  features  of  Fig.  [2]  are  the 


persistence  of  ponderomotive  scaling  for  a  1,  and  the 
even  steeper  slope  for  30  <  a  <  50.  These  features  il¬ 
lustrate  that  the  transition  between  ponderomotive  and 
phase  resonance  acceleration  is  characterized  by  a  rapid 
increase  in  energy,  rather  than  the  onset  of  a  limiting 
factor  in  the  scaling  law.  It  is  this  feature  that  makes 
xLIPA  interesting  not  only  from  the  strong  field  physics 
point  of  view,  but  also  from  the  point  of  view  of  free 
space  acceleration  of  electrons. 

Finally,  RR  has  only  a  small  effect  on  the  maximum 
photoelectron  energy,  even  for  highly  speculative  values 
of  a.  However,  the  effect  on  the  momentum  distribution 
is  interesting,  as  is  elaborated  upon  below.  RR  effects 
are  limited  to  the  multi-exawatt  regime,  and  should  play 
no  role  in  near-term  experiments. 


III.  ANALYSIS  OF  PHASE  RESONANCE 

In  this  section,  useful  estimates  are  derived  for  the 
maximum  energy  obtained  by  a  free  electron  that  is 
abruptly  introduced  into  an  extreme  laser  field.  The 
analysis  is  based  on  the  expectation  that  when  an  inner 
shell  electron  is  tunnel-ionized,  it  is  accelerated  abruptly 
to  the  speed  of  light.  In  this  limit,  the  motion  is  nearly 
parallel  to  the  wavevector  of  the  radiation,  and  the  phase 
of  the  particle  in  the  radiation  field  can  be  regarded  as 
constant.  The  primary  constraint  is  that  the  interaction 
is  limited  to  regions  where  the  irradiance  is  high  and  the 
phase  velocity  is  close  to  c.  This  corresponds  to  the  two 
regions  just  outside  the  confocal  region.  That  is,  far  from 
the  confocal  region  the  irradiance  is  too  low,  but  inside 
the  confocal  region  the  phase  velocity  is  too  high. 

There  is  an  exact  solution  for  the  motion  of  a  charged 
particle  in  any  radiation  field  of  the  form  F^u(k^x^), 
where  k ^  is  the  four  dimensional  wavevector  of  the  ra¬ 
diation,  and  x n  is  the  spacetime  coordinate.  Without 
loss  of  generality,  let  k%  =  &2  =  0.  Then 

T  =  u0  —  u3  (1) 

is  invariant  [3IE],  where  is  the  four  dimensional  ve¬ 
locity  of  the  particle1.  Combining  this  with  the  identity 
=  1,  one  can  show  that  in  the  high  energy  limit 
{ui,U2}.  In  other  words,  the  particle  moves  pre¬ 
dominate  in  the  “forward”  direction,  i.e.,  parallel  to  the 
radiation  wavevector. 

Consider  the  lowest  order,  linearly  polarized,  Hermite- 
Gaussian  laser  mode,  with  x\  the  polarization  direction 
and  xs  the  propagation  direction.  The  equations  of  mo¬ 
tion  for  a  perfectly  resonant  particle  can  be  integrated 
most  conveniently  in  the  case  X2  =  0,  so  that  the  axial 
magnetic  field  vanishes.  As  will  be  shown  below,  this  is 


1  In  particular,  u  =  (7,7^1, 7^2>7fe)>  where  eft  are  the  compo¬ 

nents  of  the  three  dimensional  velocity,  and  7  =  (1  —  /32)-1/2 
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the  most  interesting  case  for  high  energy  photoelectron 
production.  In  matrix  form,  the  equations  of  motion  are 


dx 

—  —  cu 
ds 

(2) 

du 

(3) 

—  =ftu 
ds 

where  x(s)  is  the  world  line  of  the  particle,  u(s)  is  the 
four- velocity,  and  H  =  auoF.  The  parameter  s  is  the 
proper  time,  a  =  qE/mcuo ,  E  is  the  electric  field,  q  is 
the  charge  of  the  particle,  m  is  the  mass,  and  uj  is  the 
frequency  of  the  radiation.  Using  the  coordinate  system 
described  above,  the  matrix  F  is 

0  10  e  \ 

1  0  0  -! 

0  0  0  0  1  ; 

e  1  0  0  / 

Here,  e  is  the  ratio  of  axial  to  transverse  electric  field, 
which  need  not  be  small.  Define  a  phase  resonant  par¬ 
ticle  as  one  for  which  H  is  slowly  varying  on  x(s).  Such 
particles  have  x(s)  confined  to  the  intersection  of  two  re¬ 
gions,  one  being  a  neighborhood  about  a  hypersurface  of 
constant  e,  and  the  other  a  neighborhood  about  a  hyper¬ 
surface  of  constant  phase.  In  the  first  approximation,  Q 
is  constant,  and  the  solution  of  the  velocity  equation  is 

u(s)  =  e^s?i(0)  (5) 

The  matrix  exponential  ens  is  tractable,  and  can  be  easily 
computed  using  symbolic  math  software.  Let 

A  (s)  =  ens  (6) 

It  can  be  shown  that  A  is  a  Lorentz  transformation,  i.e., 
ATgA  =  g,  where  T  indicates  the  transpose,  and  g  = 
diag(l,  —  1,  —  1,  —  1).  Of  particular  interest  is  the  initial 
condition  r(0)  =  (1,0,0,0)T,  which  according  to  most 
theories  holds  for  an  electron  at  the  moment  of  ionization, 
at  least  when  the  atomic  number  satisfies  Z  <C  137.  In 
this  case, 


/  cosh  ea  —  1  +  e2  cosh  ea 
1  e  —  e  cosh  ea  +  e  sinh  ea 
?  0 

\  cosh  ea  —  1  +  e2  sinh  ea 


(7) 


forward  direction,  i.e.,  u%  ^  u\.  Assuming  a  1,  this 
requires  that  ujs  be  at  least  of  order  unity,  i.e.,  the  time 
elapsed  according  to  a  clock  moving  with  the  particle 
should  read  at  least  one  laser  period,  as  measured  by  a 
lab  frame  clock.  This  does  not  necessarily  violate  the 
assumption  that  the  particle  should  stay  in  phase,  since 
the  two  clocks  may  keep  very  different  time. 

In  order  to  estimate  the  maximum  energy  gain,  values 
for  s  and  e  are  needed.  The  distance  traversed  by  the 
particle  is 


A 

2ire3a 


(sinh 


ea  —  ea 


e2  cosh  ea ) 


(9) 


where  A  =  27tc/cj.  The  maximum  energy  gain  is  uq(s'), 
where  s'  is  the  solution  of  xs(s')  =  zr.  Here,  zr  = 
7rro/A  is  the  Rayleigh  length,  with  ro  the  radius  of  the 
beam  waist.  A  characteristic  value  for  e  is  obtained  by 
evaluating  the  field  at  the  point  (ro,0,  zr),  which  gives 
e2  =  X/A'kzr. 

A  simple  formula  for  the  maximum  energy  gain  can 
be  obtained  by  using  the  plane  wave  estimate  for  the 
momentum.  In  this  case,  the  proper  time  is  related  to  x% 
by 


(6a;3)1/3 
(ak) 2/3 


uo{x3)  =  1  +  I  (6akx3)2/3 


so  that  the  energy  can  be  put  explicitly  in  terms  of  £3: 

(ii) 

Substituting  the  Rayleigh  length  for  X3  gives. 

Wo, max  =  1  +  1  (^)4/3  (12a)2/3  (12) 

As  an  example,  a  10  PW  laser  pulse,  with  A  =  0.8  /im, 
focused  to  ro  =  5  /im,  gives  a  «  100,  and  uo^^mc2  « 
1.5  GeV.  Comparison  of  the  prediction  accounting  for 
the  axial  field  with  that  of  Eq.  (12)  is  in  Fig.  [2]  The 
axial  field  correction  curve  is  obtained  by  solving  for  s' 
numerically. 

In  the  first  approximation,  there  is  no  advantage  in 
varying  the  laser  wavelength.  To  see  this,  rewrite  Eq.  ( 12 ) 
as 

,  2/3  /  x  i/3 


^0,n 


=  1+2(+p!)  (^t 


me 0 


(13) 


where  a(s)  =  acus.  This  gives  T  =  coshecr  —  sinhecr,  so 
that  the  invariance  of  T  is  broken.  In  the  plane  wave 
limit  (e  0)  the  particle  momentum  is 


u(s) 


1  +  cr2/ 2 
a 


0 

a2/ 2 


(8) 


and  the  invariance  of  T  is  restored,  i.e.,  T  =  1.  Moreover, 
when  a  1,  the  momentum  is  predominantly  in  the 


where  P  is  the  laser  power  and  re  is  the  classical  electron 
radius.  The  wavelength  appears  only  in  the  combination 
7*0/ A,  which  is  fixed  by  the  focusing  geometry. 

The  essential  element  in  the  forgoing  analysis  is  the  as¬ 
sumption  that  upon  ionization  into  an  extreme  field,  an 
electron  can  be  accelerated  to  nearly  the  speed  of  light 
in  a  fraction  of  an  optical  cycle.  This  requires  that  the 
ionization  potential  be  large  enough  so  that  the  electron 
is  held  in  position  by  its  parent  ion  until  it  is  exposed  to 
ultra-relativisitic  intensity,  but  not  so  large  that  ioniza¬ 
tion  becomes  highly  improbable. 
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IV.  TWO-STEP  IONIZATION  MODEL 

An  elegant,  but  incomplete,  calculation,  of  the  mo¬ 
mentum  distribution  of  relativistic  tunnel-ionized  elec¬ 
trons,  can  be  carried  out  in  the  strong  field  approxima¬ 
tion  (SFA)  p~6j.  Recently,  the  SFA  was  extended  to  ac¬ 
count  for  the  Coulomb  field  of  the  residual  ion  by  means 
of  a  quasi-classical  correction  factor  [6] .  An  alternative 
approach  is  the  imaginary  time  method  (ITM)  [8].  The 
SFA  is  perhaps  more  suitable  for  the  purposes  of  the 
present  paper,  due  to  the  fact  that  it  utilizes  the  S- 
matrix,  which  gives  directly  the  primary  object  of  inter¬ 
est,  the  momentum  distribution  of  the  ionized  electron. 
For  this  reason,  the  following  discussion  emphasizes  the 
SFA  over  the  ITM. 

A  major  claim  of  this  paper  is  that  the  assumptions 
typically  utilized  in  practical  SFA  analyses  fail  in  the 
extreme  light  regime.  Foremost  among  these,  are  the  re¬ 
strictions  on  the  form  of  the  laser  field,  which  is  invariably 
taken  to  have  perfectly  planar  phase  fronts,  and  slowly 
varying  envelope.  In  order  to  remove  these  assumptions, 
the  two-step  model  suggested  by  Corkum  [2]  is  employed. 
In  the  first  step,  a  free  electron  is  spawned  at  a  random 
time,  chosen  so  that  the  probability  of  ionization  during 
any  time  interval  is  consistent  with  some  rate  law  (see 
appendix).  In  the  second  step,  the  free  electron  evolves 
according  to  the  classical  laws  of  motion  until  it  leaves 
the  interaction  region.  Carrying  out  a  large  number  of 
trials  gives  what  may  be  termed  a  classical  S-matrix,  Srp , 
where  |  ^rp  | 2  is  the  probability  that  an  ion  with  spatial  co¬ 
ordinate  r  will  produce  a  photoelectron  with  momentum 
p.  This  dependence  on  the  ion  coordinates,  absent  in  the 
SFA,  comes  about  because  of  finite  spot-size  effects. 

The  two-step  model  requires  a  rate  law,  which  in  this 
work  is  the  Coulomb  corrected,  dressed,  SFA  of  Klaiber 
et  al.  [6] .  This  creates  the  difficulty  that  a  subset  of  the 
information  contained  within  the  SFA  (the  rate  law)  is 
used  to  obtain  a  result  that  differs  from  the  total  infor¬ 
mation  contained  within  the  SFA  (the  momentum  distri¬ 
bution).  In  order  for  this  to  make  sense,  it  is  required 
that  the  results  from  the  two-step  model  should  not  de¬ 
pend  strongly  on  which  rate  law  is  used.  To  this  end,  we 
verified  that  the  momentum  distributions  derived  from 
the  two-step  model  are  nearly  indistinguishable  whether 
the  ITM,  the  dressed  SFA,  or  the  undressed  SFA,  are 
employed  in  the  rate  law. 

Confirmation  that  the  two-step  model  agrees  with  the 
SFA,  in  the  limit  of  a  plane  wave,  is  given  in  Fig.  [3j  The 
peak  normalized  vector  potential  of  the  laser  is  a  =  36, 
and  the  pulse  width  is  20  fs.  The  ion  under  consideration 
is  hydrogen- like  argon,  Z  =  18.  Various  projections  of 
the  differential  ionization  rate  are  plotted  using  both  the 
SFA  and  the  two-step  model.  In  computing  a  differential 
ionization  rate  from  the  two-step  model  output,  the  S- 
matrix  probabilities  have  to  be  divided  by  a  suitable  time 
interval.  The  time  interval  used  to  generate  Fig.  [3]  is 
the  width  of  the  ionization  current  pulse  envelope.  This 
is  determined  by  monitoring  the  number  of  ionization 
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FIG.  3:  Benchmarking  the  two-step  model  in  the  plane  wave 
limit,  with  a  =  36,  Z  =  18,  and  r  =  20  fs,  by  comparing  the 
differential  ionization  rate  with  the  SFA  prediction.  Panel  (a) 
is  the  analytical  SFA  formula,  (b)  is  the  numerical  two-step 
model,  and  (c)  projects  out  U3  to  allow  for  better  quantitative 
comparison.  In  (c),  the  blue  curve  is  the  SFA,  while  the  green 
curve  is  the  two-step  model. 


events  per  unit  time,  integrated  over  all  trials. 

The  parabolic  shape  of  the  distribution  in  the  ui~us 
plane,  which  follows  from  the  invariance  of  T,  is  indis¬ 
tinguishable  from  one  case  to  the  other.  The  relative 
density  distribution  along  the  parabola  is  also  very  simi¬ 
lar.  Comparison  of  the  absolute  scale  in  the  ui-u%  plane 
is  not  meaningful,  because  strictly,  the  two-step  model 
distribution  contains  a  factor  S(T  —  1),  whereas  the  SFA 
allows  for  small  deviations  about  T  =  1.  The  value  of 
dw / duidus  is  rendered  finite  in  Fig.  [3^b)  by  applying  a 
smoothing  filter.  In  order  to  make  a  comparison  where 
the  absolute  scale  has  meaning,  Fig.  [3jc)  shows  the  dif¬ 
ferential  ionization  rate  projected  onto  the  u\  axis. 

Fig.  [^establishes  that  the  two-step  model  is  capable  of 
generating  a  momentum  distribution  that  closely  approx- 
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FIG.  4:  Representative  orbits  of  photoelectrons  for  various 
parent  ions.  Electrons  from  the  M  or  L  shell  exhibit  pondero- 
motive  acceleration,  while  electrons  from  the  K  shell  exhibit 
phase  resonance.  The  fractional  power  scale  on  the  horizon¬ 
tal  axis  gives  the  effect  of  a  log  scale,  while  allowing  for  a 
transition  through  the  origin. 


imates  the  SFA  in  the  plane  wave  limit.  An  important 
characteristic  of  the  plane  wave  case  is  that  the  initial 
position  of  the  parent  ion  does  not  affect  the  final  mo¬ 
mentum  in  any  way.  The  two-step  model  allows  finite 
spot-size  illumination  to  be  considered,  and  the  effect  of 
the  ion  coordinates  to  be  determined. 
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V.  MOMENTUM  DISTRIBUTIONS 

Consider  a  system  comprised  of  many  ions  distributed 
in  space.  Assuming  that  correlations  among  photoelec¬ 
trons  are  negligible,  the  total  momentum  distribution  is 

/( p)  =  ]T|Srip|2  (14) 

i 

where  i  indexes  each  ion.  By  loading  a  large  number  of 
ions  into  any  given  simulation,  an  ensemble  of  ion  co¬ 
ordinates  and  ionization  times  is  created,  and  /( p)  is 
generated  by  binning  the  whole  set  of  photoelectrons  in 
momentum  space. 

Consider  laser  parameters  a  =  100,  A  =  0.8  /im,  ro  =  5 
fi m,  and  r  =  20  fs,  which  are  chosen  to  be  achievable 
by  near-term  10  PW  laser  systems.  Representative  pho¬ 
toelectron  orbits  from  the  two-step  model  are  shown  in 
Fig-@  for  various  parent  ions.  The  optimum  energy  gain 
is  obtained  for  Ar17+.  The  lower  charge  states  of  argon 
are  ionized  too  early  in  the  foot  of  the  pulse,  and  are 
ponderomotively  ejected  before  phase  resonance  can  be 
reached.  The  case  of  Ti21+  is  similar  to  that  of  Ar17+, 
while  that  of  Fe25+  gives  more  limited  energy  gain.  The 
limitation  in  the  case  of  Fe25+  comes  about  because  the 
photoelectrons  tend  to  be  spawned  too  near  the  origin, 
which  is  not  the  ideal  starting  location. 

The  momentum  distribution,  /( p),  of  photoelectrons 
extracted  from  a  uniform  distribution  of  Ar17+  ions,  is 
shown  in  Fig.  [5j  Comparison  of  Fig.  |5ja)  with  Figs.  la) 


FIG.  5:  Projections  of  the  photoelectron  momentum  distribu¬ 
tion,  /( p),  for  an  ensemble  of  initial  conditions  corresponding 
to  uniformly  distributed  Ar17+  ions. 

or  (b)  reveals  that  finite  spot  size  effects  lead  to  a  pro¬ 
found  change  in  the  photoelectron  spectrum.  This  can  be 
easily  understood  in  terms  of  the  fact  that  in  the  plane 
wave  case,  there  is  no  possibility  of  any  electron  stay¬ 
ing  confined  to  a  small  region  of  phase  throughout  the 
interaction.  The  interaction  stops  only  after  all  phases 
of  the  plane  wave  pass  by.  When  the  spot  size  is  finite, 
the  interaction  is  terminated  when  the  electron  exits  the 
confocal  region.  In  the  latter  case,  phase  resonance  is 
possible,  and  higher  energies  are  obtained. 

Projections  of  the  classical  S-matrix  for  Ar17+,  corre¬ 
lating  final  energy  with  initial  coordinate,  are  shown  in 
Fig.©  In  any  given  projection,  there  is  a  large  low  energy 
population,  and  only  a  small  number  of  high  energy  parti¬ 
cles.  The  highest  energy  particles  originate  roughly  from 
the  coordinates  x\  ~  ±t*o,  x 2  ~  0,  and  X3  ~  —2 zr.  This 
suggests  that  a  high  quality,  high  energy  beam,  might  be 
obtained  by  localizing  the  parent  ions  to  a  small  neigh¬ 
borhood  about  one  or  both  of  these  two  points.  If  a  gas 
target  is  used,  some  advantage  might  be  gained  by  fo¬ 
cusing  into  a  plasma  channel  or  lens  mm,  so  that  the 
highest  ion  density  is  off-axis.  Alternatively,  one  might 
consider  suspending  titanium  nanoparticles  at  the  appro¬ 
priate  points  in  the  laser  focus.  Highly  accurate  position¬ 
ing  and  alignment  would  be  needed  in  either  case,  and 
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FIG.  6:  Classical  S-matrix  projections  for  Ar17+.  The  cor¬ 
relation  of  final  energy  with  initial  polarization  coordinate  is 
in  (a),  with  initial  cross-polarization  coordinate  is  in  (b),  and 
with  initial  axial  coordinate  is  in  (c). 


pre-pulses  would  have  to  be  strongly  extinguished. 

Leaving  these  practical  matters  in  abeyance,  consider 
the  photoelectron  distribution  from  a  single  Ar17+  ion  po¬ 
sitioned  near  the  optimal  coordinate.  This  is  computed 
using  106  trials  of  the  two-step  model,  with  ion  coor¬ 
dinate  x  =  (t*o,0,  —  2zr).  The  results  are  illustrated  in 
Fig.  [7|  The  ui-us  projection,  shown  in  Fig.  [TJa),  shows 
a  clear  preference  for  positive  u\  values,  evidently  due 
to  the  choice  of  positive  x\.  The  fine  structure  exhibits 
discrete  filamentary  structures,  including  both  open  and 
closed  loops.  This  structure  can  be  understood  by  in¬ 
specting  Fig.  [7]^b),  which  shows  the  S'-matrix  projected 
into  the  plane  of  ionization  phase  and  energy.  The  ioniza¬ 
tion  phase  is  defined  as  ip  =  u(tc  —  x$/c),  where  tc  is  the 
time  of  ionization  (see  appendix).  The  transverse  electric 
field  component,  evaluated  at  the  point  x,  is  plotted  as 
a  function  of  ip  in  Fig.  (7^c).  Because  of  finite  spot-size 
and  pulse- width  effects,  the  peaks  occur  at  non-integral 
values  of  p/ir.  Cross-referencing  Figs.  13 a),  (b),  and  (c) 
makes  it  apparent  that  the  closed  loops  in  (a)  are  gen¬ 
erated  by  particles  ionized  near  the  peaks  of  the  electric 
field.  The  open  loops  are  generated  either  just  before 
the  electric  field  reaches  its  peak  negative  value,  or  just 
after  it  reaches  its  peak  positive  value.  The  correspond¬ 
ing  electrons  are  situated  to  stay  in  phase  resonance  for 
the  longest  time,  and  therefore  reach  the  highest  ener¬ 
gies.  Finally,  the  energy  spectrum  within  a  cone  angle 
9  =  0.01  radians  is  shown  in  Fig.  (7^d).  The  three  narrow 
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FIG.  7:  Momentum  distribution  from  a  single,  well  posi¬ 
tioned,  Ar17+  ion.  The  u\-us  momentum  projection  is  in 
(a),  and  the  ip-uo  ^-matrix  projection  is  in  (b).  The  value  of 
Ei,  evaluated  at  the  position  of  the  ion,  is  shown  in  (c)  as  a 
function  of  p.  The  energy  spectrum  for  photoelectrons  inside 
a  10  milliradian  cone  angle  is  in  (d). 


peaks  around  1.5  GeV  are  associated  with  the  resonant 
electrons  that  originate  in  positive  electric  fields.  The 
broad,  higher  energy  tail  is  due  to  the  resonant  electrons 
that  originate  in  negative  electric  fields. 
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FIG.  8:  Momentum  distribution  from  a  single,  well  posi¬ 
tioned,  Au78+  ion,  exposed  to  an  irradiance  of  2  x  1026 
W/cm2.  The  case  neglecting  RR  is  in  (a),  and  the  case  with 
the  exact  Landau  and  Lifshitz  formula  for  RR  is  in  (b). 


VI.  RADIATION  REACTION  EFFECTS 


to  the  SFA,  is  sufficient  to  fully  strip  gold  (Z  =  79).  Con¬ 
sider  placement  of  an  Au78+  ion  at  the  same  location  in 
the  confocal  region  considered  above,  x  =  (ro,0,  —2 zr). 
The  u\-uz  momentum  distributions  with  and  without  RR 
are  displayed  in  Fig.  [8]  In  this  extreme  case,  the  charge 
is  peaked  at  the  highest  energies,  even  without  any  angu¬ 
lar  selection.  Either  with  or  without  RR,  there  are  two 
main  peaks  in  the  energy  spectrum.  The  peaks  in  the 
spectrum  with  RR  are  lower  energy,  but  narrower.  The 
angular  spread  in  the  case  with  RR  is  also  reduced.  The 
overall  shape  in  the  two  cases  is  obviously  quite  different. 
Developing  a  full  understanding  of  these  effects  promises 
to  be  a  rich  area  of  research. 


VII.  CONCLUSIONS 


Multi-petawatt  lasers  access  a  new  regime  of  strong 
field  physics  and  free  space  acceleration,  which  appears 
when  the  normalized  vector  potential  a  >  100.  The  elec¬ 
tric  field  in  the  confocal  region  is  sufficient  to  fully  strip 
moderately  heavy  atoms  such  as  argon  or  titanium  via 
tunneling  ionization.  The  K-shell  electrons  are  released 
into  such  a  large  electric  field,  they  are  accelerated  to 
the  speed  of  light  in  a  fraction  of  an  optical  period,  and 
are  deflected  into  the  laser  propagation  direction  by  the 
magnetic  field.  Some  of  these  electrons  are  accelerated  to 
gigavolt  energies  due  to  phase  resonance,  i.e.,  they  stay 
near  the  same  optical  phase  throughout  a  substantial  por¬ 
tion  of  the  confocal  region.  The  resulting  photoelectron 
momentum  distributions  have  unique  features,  which  de¬ 
part  significantly  from  the  usual  parabolic  form  predicted 
by  the  SFA  and  other  analyses.  For  atoms  loaded  into 
a  small  spatial  region,  the  momentum  distribution  has 
discrete  features  that  can  be  related  to  sub-cycle  bursts 
of  ionization  current. 


A  fundamental  problem  in  the  strong  field  physics  of 
atomic  systems,  which  seems  to  have  been  left  unexplored 
up  until  now,  is  the  effect  of  radiation  reaction  (RR)  on 
the  spectrum  of  tunnel  ionized  electrons.  RR  is  the  force 
acting  on  an  electron  due  to  its  own  fields.  This  force 
is  normally  negligible,  only  becoming  significant  for  very 
large  fields  or  electron  energies.  With  laser  fields  reach¬ 
ing  ever  more  stupendous  values,  the  problem  of  RR  is 
attracting  renewed  interest  [MSI.  The  fields  needed 
to  bring  the  photoelectron  dynamics  into  the  RR  dom¬ 
inated  regime  are  far  in  excess  of  what  is  likely  to  be 
achievable  in  the  near  future.  This  can  be  understood 
in  terms  of  the  fact  that  high  energy  photoelectrons  al¬ 
ways  co-propagate  with  the  optical  wave,  whereas  RR  is 
encouraged  by  counter-propagating  geometries. 

The  numerical  methods  used  here  incorporate  the  ex¬ 
act  Landau  and  Lifshitz  formula  for  RR  (see  appendix). 
Significant  effects  of  RR  on  the  photoelectron  spectrum, 
given  the  focal  geometry  considered  throughout  this  pa¬ 
per,  begin  to  appear  for  a  =  5000,  corresponding  to  about 
25  exawatts  of  laser  power.  Such  a  laser  field,  according 
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Appendix  A:  Extreme  Field  Particle  Pusher 

Ordinary  particle  pushers  are  ineffective  in  extreme 
fields.  In  Ref.  [25],  a  covariant  pusher  with  an  adaptive 
time  step  is  used  to  maintain  accuracy.  Here,  the  co¬ 
variant  form  is  retained,  but  an  important  improvement 
is  made  which  eliminates  the  need  for  explicit  time  step 
adjustment,  and  increases  computational  efficiency  by  or- 


ders  of  magnitude.  The  treatment  of  radiation  reaction 
is  identical  to  the  one  in  03  • 

The  exact  operator  of  momentum  evolution  is  a 
Lorentz  transformation.  Typically  one  approximates  this 
as  some  sequence  of  boosts  and  rotations  [25ti27j.  Due 
to  the  fact  that  these  are  not  commuting  operators,  an 
error  is  introduced,  which  is  an  increasing  function  of  the 
field  strength.  To  avoid  this  problem,  a  compact  form  of 
the  complete  Lorentz  transformation  is  needed. 

The  pusher  used  here  advances  particles  in  proper 
time,  and  takes  a  manifestly  covariant  form.  One  use¬ 
ful  feature  of  a  proper  time  advance  is  that  a  constant 
proper  time  step  leads  to  a  constant  phase  step  for  any 
particle  in  a  plane  wave,  no  matter  its  energy  or  direction. 
In  fact,  requiring  that  the  phase  step  satisfy  Acp  27 r 
leads  to 


A,  «  -  (15) 

where  T  is  a  constant  of  the  motion.  For  a  particle  ini¬ 
tially  at  rest,  T  =  1,  and  the  proper  time  step  appro¬ 
priate  for  any  particle  is  the  same  as  the  lab  frame  time 
step  appropriate  for  a  non-relativistic  particle. 


The  condition  (15)  is  sufficient  provided  the  momen¬ 
tum  advance  is  exact  in  a  uniform,  constant  field.  This 
requirement  is  satisfied  by  the  expression  u(s  +  As/2)  = 
A(s,  As)u(s  —  As/2),  where  s  is  the  proper  time,  As  is 
the  time  step,  and 


A(s,  As)  =  exp  (qF(s)As/mc)  (16) 

Here,  F(s)  is  the  antisymmetric  field  tensor  evaluated  on 
the  world  line  x(s).  To  avoid  the  appearance  of  oner¬ 
ously  complicated  expressions,  it  is  convenient  to  first 
rotate  the  coordinate  system  so  that  there  is  only  one 
component  of  B,  apply  the  matrix  exponential,  and  fi¬ 
nally  restore  the  coordinate  system.  Specifically, 

A (s,  As)  =  1  1(B)Ao(s,  As)T(B)  (17) 

where  T (B)  aligns  the  magnetic  field  with  one  of  the  basis 
vectors.  For  definiteness  let  this  be  the  basis  vector  in 
the  1-direction,  ei.  Then,  T(B)  is  any  rotation  matrix 
satisfying  T(B)  •  B  =  |B|ei.  The  matrix  exponential 
is  conveniently  split  into  symmetric  and  antisymmetric 
parts,  such  that  Aq  =  A^  +  A  a,  where 
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Here,  dotted  entries  are  used  to  emphasize  that  the  cor¬ 
responding  matrix  elements  can  be  obtained  from  sym¬ 
metry.  The  several  variables  appearing  in  the  matrix 
exponential  are  defined  in  Table  [TJ  In  programming  the 
matrix,  care  must  be  taken  to  prevent  floating  point  ex¬ 
ceptions  in  regions  where  a  or  a±  vanish,  such  as  in  field 
free  regions,  or  ideal  plane  waves.  Imposing  a  miniscule 
uniform  field  is  one  simple  solution.  Also,  some  advan¬ 
tage  might  be  gained  by  updating  one  of  the  components 
of  u  via  the  identity  uT gu  =  1,  rather  than  by  using  the 
corresponding  row  of  A0.  In  this  work  uq  is  so  updated. 

In  order  to  emphasize  the  importance  of  the  covariant 
pusher  just  described,  its  accuracy  is  compared  with  that 
of  the  Boris  pusher  m  in  Fig.  [9j  The  simulation  initial¬ 
izes  a  single  electron  into  a  monochromatic  plane  wave 
such  that  T  =  1  and  u\  =  —qAi/mc2.  The  error  mea¬ 
sure  1  —  T,  which  should  be  zero  at  all  times,  is  plotted 
in  Fig.  [9|a).  The  accuracy  of  the  Boris  pusher  suffers 


with  increasing  a,  while  the  covariant  pusher  is  accurate 
for  any  a.  Phase  space  orbits  for  a  =  10  are  plotted  in 
Fig.  gb).  The  covariant  pusher  correctly  recovers  the 
spatial  period,  which  is  known  analytically  to  be  a2 A/4. 
In  contrast,  the  Boris  pusher  gives  a  spuriously  long  spa¬ 
tial  period. 

It  is  important  to  note  that  the  covariant  pusher  gives 
accurate  results  much  faster  than  the  Boris  pusher.  The 
time  step  used  in  Fig.  [9]  is  27t/267cj.  In  the  case  of  the 
Boris  pusher,  the  time  step  is  measured  in  the  lab  frame, 
whereas  in  the  case  of  the  covariant  pusher,  it  is  measured 
in  the  frame  of  the  electron.  The  time  elapsed  in  the  lab 
frame  is  always  longer  than  the  time  elapsed  according  to 
the  electron’s  own  clock.  As  a  result,  the  covariant  pusher 
always  requires  fewer  steps  than  the  Boris  pusher.  This 
is  related  to  the  fact  that  the  phase  step  in  the  covari¬ 
ant  pusher  is  constant,  while  the  phase  step  in  the  Boris 
pusher  decreases  with  increasing  forward  momentum.  As 
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TABLE  I:  Variables  Appearing  in  Aq 


Variable 

Definition 

(TE)  •  e; 

b 

(TB)  •  ei 

e2 

el  +  e2  +  e3 

a2 

4 b2el  +  (e2  -  b2)2 

a± 

a  d=  (e2  —  b2) 

X+ 

cosh(o+^As/v/2nrc) 

X- 

cos  (a  _  q  A  s/ \f2m  c) 

CG 

sinh(a+^As  /  y/2mc)  / 

<j- 

sin(o_  q  As  /  y/2  me)  /a_ 

(X+  ±X-)/2 

E± 

(cr+  ±cr_)/2 

F 

b2  +  e2 

G 

—b2  —  e\  +  el  +  el 

H 

b2  +  e2  —  e  2  +  e2 

J 

b2  +  e\  +  e2  —  e2 

Boris  Pusher 
o  Covariant  Pusher 


(a) 


(b) 


an  example,  in  Fig.  |9jb),  the  covariant  pusher  curve  is 
generated  in  267  steps,  while  the  Boris  pusher  curve  is 
generated  in  6800  steps. 

It  is  natural  to  ask  whether  using  an  adaptive  time  step 
is  an  alternative  to  using  the  unsplit  matrix  exponential 
A.  This  approach  is  elaborated  upon  in  [25j,  where  the 
pusher  is  in  covariant  form,  but  retains  the  conventional 
time  step  splitting.  This  turns  out  to  be  far  less  efficient 
than  using  the  unsplit  A. 

Another  advantage  of  expressing  the  particle  pusher 
in  covariant  fashion  is  that  the  Lorentz-Abraham-Dirac 
formula  for  the  radiation  reaction  force  takes  the  simple 
form 


R  = 


2  q2 

3  me 


d2u 
ds 2 


T 

uu  g 


d2u 
ds 2 


(20) 


The  Landau  and  Lifshitz  (LL)  formula  [28j  is  derived  by 
substituting  for  d2u/ds 2  the  value  obtained  in  the  ab¬ 
sence  of  radiation  reaction.  In  three  dimensional  nota¬ 
tion,  the  LL  formula  is  extremely  unwieldy,  and  even  in 
four  dimensional  form,  it  appears  to  require  expensive 
evaulations  of  all  spacetime  derivatives  of  the  field  ten¬ 
sor.  When  the  covariant  pusher  described  above  is  used, 
a  simple  and  elegant  alternative  becomes  readily  avail¬ 
able.  Namely,  by  splitting  each  step  into  two  half-steps, 
d2u/ds 2  can  be  evaluated  by  direct  finite  differencing. 
That  is,  during  each  step  generate 

u(1)  =  a(*+T’t)“(0)  (21) 


i.(2) 


.  3 As  As\  (1) 

A(s+  — ,—  U(1) 


(22) 


This  requires  field  evaluations  at  only  two  spacetime 
points.  Now,  in  the  absence  of  radiation  reaction,  = 
u(s),  u W  =  u(s  +  As/2),  and  u ^  =  u(s  +  As).  The  LL 
reaction  force  is  therefore  given  by  the  matrix  equation 

R  =  — —  (S2u  —  uuTgS2u)  (23) 

3 me  v  7 

where 


S2u  = 


u —  2 +  u 
As2/4 


(24) 


is  the  finite  difference  form  of  dPu/ds2.  The  updated 
four- velocity,  including  the  LL  reaction  force,  is 

u(s  +  As)  =  u ^  +  RAs  (25) 


Here,  a  simple  Euler  advance  is  justified  by  the  expecta¬ 
tion  that  the  reaction  force  is  small. 


FIG.  9:  Errors  due  to  the  Boris  pusher  and  the  covariant 
pusher  introduced  in  this  paper.  The  error  measure,  1  —  T,  is 
plotted  as  a  function  of  a  in  (a) ,  and  the  orbits  corresponding 
to  the  case  a  —  10  are  plotted  in  (b). 


Appendix  B:  Elaboration  on  the  Two-Step  Model 

The  quantum  mechanical  S'-matrix  is  at  the  foundation 
of  the  SFA.  It  has  a  simple  interpretation  as  the  matrix  of 
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amplitudes  of  all  possible  transitions  between  in  and  out 
states.  The  ^S-matrix  viewpoint  holds  that  the  electron 
dynamics  during  the  interaction  are  not  observable.  The 
S'-matrix  expresses  this  fundamental  ignorance  by  giving 
only  the  transition  amplitude  between  states  in  the  dis¬ 
tant  past  or  future,  when  the  interaction  is  negligible.  In 
contrast,  the  two-step  model  allows  for  a  full  dynamical 
description  of  the  electron  motion  after  a  certain  critical 
point  in  time,  which  is  loosely  identified  with  the  “time 
of  ionization.”  This  raises  a  problem  of  interpretation 
due  to  the  fact  that  time  is  not  a  quantum  mechanical 
observable. 

Before  addressing  the  resolution  of  this  problem,  it  is 
useful  to  note  that  many  quantum  mechanical  ^S-matrix 
calculations  do  in  fact  introduce  something  like  a  critical 
time  in  the  course  of  making  mathematical  approxima¬ 
tions.  In  particular,  the  total  S'-matrix  amplitude  is  ob¬ 
tained  from  the  saddle  point  method,  which  introduces  a 
summation  over  saddle  points,  each  of  which  corresponds 
to  a  particular  phase  of  the  applied  field.  There  are  two 
saddle  points  per  optical  cycle.  One  is  tempted  to  in¬ 
terpret  this  procedure  as  a  summation  over  ionization 
events,  which  are  localized  in  phase  near  the  peaks  of  the 
field. 

In  order  to  interpret  the  time  of  ionization,  Bohmian 
mechanics  is  used.  Bohmian  mechanics  is  equivalent  to 
quantum  mechanics  in  all  its  statistical  predictions,  but 
treats  the  electron  coordinate  as  a  hidden  variable  [29]. 
The  electron  responds  to  electromagnetic  forces  classi¬ 
cally,  but  also  responds  to  a  quantum  force  that  depends 
on  the  wavefunction.  In  this  context,  the  time  of  ioniza¬ 
tion  is  easy  to  define.  Let  tc  be  the  earliest  time  such  that 
for  all  t  >  £c,  the  classical  force  dominates  the  quantum 
force.  Then  tc  is  the  critical  time  after  which  the  electron 
behaves  classically.  This  is  as  useful  a  notion  of  “time  of 
ionization”  as  one  is  likely  to  discover. 

In  the  Bohmian  picture,  the  coordinates  of  an  elec¬ 
tron  in  a  stationary  state  can  take  a  continuous  range  of 
definite  values.  The  statistical  distribution  of  these  ini¬ 
tial  coordinates  determines  all  the  subsequent  statistics, 
such  as  the  distribution  of  ionization  times.  In  order  to 
proceed,  it  is  not  necessary  to  solve  the  whole  Bohmian 
mechanics  problem.  Instead,  a  single  hidden  variable, 
iL,  is  introduced  into  the  initial  electron  distribution.  It 
turns  out  that  if  the  electron  distribution  contains  a  fac¬ 
tor  e~H ,  and  if  tc  is  defined  by 


trajectories  would  be  needed.  An  interesting  task  for  fu¬ 
ture  research  would  be  to  utilize  Bohmian  trajectories, 
computed  numerically  [30] .  as  inputs  into  the  two-step 
model.  In  this  work,  appeal  is  made  to  the  ITM,  which 
suggests  that  u(tc)  ~  (1,0,  0,  0)T.  The  spatial  coordinate 
of  the  ionized  electron  is  assumed  to  coincide  with  that 
of  the  parent  ion. 

In  carrying  out  the  classical  particle  tracking,  expres¬ 
sions  for  the  laser  field  components  are  needed  as  inputs 
into  the  particle  pusher  algorithm.  Let  ijj  be  the  lowest 
order  Hermite-Gaussian  solution  of  the  scalar  paraxial 
wave  equation  m-  In  order  to  impose  a  temporal  en¬ 
velope,  let  A i  =  'ipex.pl—  (z  —  vt)2 /c2t 2],  where  v  is  the 
velocity  of  the  envelope,  taken  to  be  c  in  this  work.  Then, 
to  lowest  order  in  A/ro,  A3  is  found  from  V  •  A  =  0.  As 
usual,  A2  is  neglected.  The  field  tensor  is  obtained  by 
differentiating  AM,  also  to  lowest  order  in  A/ro- 

Finally,  we  briefly  comment  on  the  electron  spin.  The 
spin  of  the  electron  is  accounted  for  in  the  SFA  rate  law, 
which  gives  the  spin  averaged  ionization  rate  [6|.  Spin 
effects  in  the  context  of  the  SFA  are  explicitly  consid¬ 
ered  in  Ref.  [32].  In  the  two-step  model,  spin  effects  are 
neglected  for  t  >  tc.  It  is  possible  to  isolate  the  spin 
terms  in  the  Dirac  equation  [33],  and  form  the  ratio  of 
spin  terms  to  other  terms.  In  a  laser  field,  this  ratio  is 
of  order  huj/amc2 ,  which  is  miniscule  for  the  parameters 
considered  in  this  paper. 


then  w(t)  is  the  ionization  rate  [25].  Therefore,  H  is  a 
hidden  variable  that  recovers  any  quantum  mechanical 
ionization  rate,  w(t),  that  one  chooses  to  use.  Moreover, 
the  time  tc  is  the  time  of  ionization,  after  which  clas¬ 
sical  particle  tracking  is  appropriate.  What  is  missing 
is  a  clear  specification  of  the  electron  coordinates  at  the 
moment  of  ionization,  x(tc)  and  u(tc).  In  order  to  rig¬ 
orously  specify  these,  knowledge  of  all  possible  Bohmian 
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